/*
 * CorrelationCoefficient.cpp
 *
 *  Created on: Nov 17, 2010
 *      Author: mgibson
 */

#include "Coefficienter.h"
#include "DataReader.h"
#include <iostream>
#include <fstream>
#include <cmath>
#include <stdlib.h>
#include "globals.h"

using namespace std;

Coefficienter::Coefficienter() {
	numberOfCrystals = 2760;
}

Coefficienter::~Coefficienter() {
	;
}


float** Coefficienter::FindCoefficients(int** packet){
	float** coefficients = new float*[numberOfCrystals];
	for(int i = 0 ; i<numberOfCrystals ; i++){
		coefficients[i] = NULL;
	}

	for(int i = 0 ; i<numberOfCrystals ; i++){
		if(packet[i]!=NULL){
			coefficients[i] = new float[48+referenceSize-1];
			for(int j = 0 ; j<48+referenceSize-1 ; j++){
				coefficients[i][j] = 0;
			}
			double packetMax = 0;
			for(int j = 0 ; j<48 ; j++){
				if(packetMax<packet[i][j]-packet[i][50])
					packetMax = packet[i][j]-packet[i][50];
			}

			double ref[referenceSize];
			for(int j = 0 ; j<referenceSize ; j++)
				ref[j] = packetMax*reference[j];

//			if(packetMax < 10) continue;

			double dataSum = 0;
			double data2Sum = 0;
			double refSum = 0;
			double ref2Sum = 0;

			for(int j = 0 ; j<48+referenceSize/2 ; j++){
				double dataRefSum = 0;
				int overlap;
				if(j<referenceSize){
					overlap = j+1;
					refSum += ref[referenceSize-overlap];
					ref2Sum += ref[referenceSize-overlap]*ref[referenceSize-overlap];
					dataSum += packet[i][j]-packet[i][50];
					data2Sum += (packet[i][j]-packet[i][50])*(packet[i][j]-packet[i][50]);

					for(int k = 0 ; k<overlap ; k++)
						dataRefSum += (packet[i][k]-packet[i][50])*ref[referenceSize-overlap+k];
				}
				else if(j>47){
					overlap = referenceSize-j+48-1;
					refSum -= ref[overlap];
					ref2Sum -= ref[overlap]*ref[overlap];
					dataSum -= packet[i][j-referenceSize]-packet[i][50];
					data2Sum -= (packet[i][j-referenceSize]-packet[i][50])*(packet[i][j-referenceSize]-packet[i][50]);

					for(int k = 0 ; k<overlap ; k++)
						dataRefSum += (packet[i][48-overlap+k]-packet[i][50])*ref[k];
				}
				else{
					overlap = referenceSize;
					dataSum += packet[i][j]-packet[i][50];
					data2Sum += (packet[i][j]-packet[i][50])*(packet[i][j]-packet[i][50]);
					dataSum -= packet[i][j-referenceSize]-packet[i][50];
					data2Sum -= (packet[i][j-referenceSize]-packet[i][50])*(packet[i][j-referenceSize]-packet[i][50]);

					for(int k = 0 ; k<referenceSize ; k++){
						dataRefSum += (packet[i][j-referenceSize+k]-packet[i][50])*ref[k];
					}
				}

				if(j > referenceSize-1)
					coefficients[i][j-(referenceSize-1)/2] = (overlap*dataRefSum-dataSum*refSum)/pow(((overlap*data2Sum-dataSum*dataSum)*(overlap*ref2Sum-refSum*refSum)),.5);
			}
		}
	}

	return coefficients;
}
